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Abstract 

Lateral distributions for electrons and muons in extensive air showers measured with 
the array of the KASCADE experiment are compared to results of simulations based 
on the high-energy hadronic interaction models QGSJet and SIBYLL. It is shown, 
that the muon distributions are well described by both models. Deviations are found 
for the electromagnetic component, where both models predict a steeper lateral 
shape than observed in the data. For both models the observed lateral shapes of the 
electron component indicate a transition from a light to a more heavy composition 
of the cosmic ray spectrum above the knee. 
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1 Introduction 



The particle lateral distribution of extensive air showers (EAS) is the key 
quantity for cosmic ray ground observations, from which most shower observ- 
ablcs are derived. The interaction cascade, which is initiated by a high energy 
cosmic ray particle in the atmosphere, creates a multitude of secondary par- 
ticles, which arrive nearly at the same time but distributed over a large area 
perpendicular to the direction of the original particle. The disc of secondary 
particles may extend over several hundred meters from the shower axis, with 
maximum density in the center of the disc, which is called the shower core. 
Apart from the arrival times, the density distribution of particles within the 
shower disc contain all informations on the primary particle, which are left 
after it has undergone a millionfold multiplication process by the atmosphere. 
However, it is this multiphcation process, that foremost offers the chance to 
observe cosmic rays in the ultra and very high energy region at all: Due to 
their low flux measurements at ground, carried out with large arrays of indi- 
vidual detectors which take samples of the shower disc at several locations, 
are still the only possible way to study these high-energy cosmic particles [1]. 

The lateral distributions of electrons and muons in EAS not only contains 
information on the nature of the primary cosmic ray particle, which is related 
to astrophysical questions, but it also carries information relevant to particle 
physics. While the electromagnetic interactions are thought to be well under- 
stood, this is not true for the high energy hadronic interactions. The energy 
range and the kinematical region in which the first hadronic interactions of 
the shower development occur are far beyond the accessible realm of todays 
accelerator experiments. Uncertainties in the description of hadronic interac- 
tions therefore imply uncertainties in the prediction of the shape of the lateral 
distributions [2]. 

A parameter, commonly used to describe the form of the lateral density dis- 
tribution, is the lateral form parameter in the Nishimura-Kamata-Greisen 
(NKG)-function [3,4,5], usually called age. The name expresses the relation 
between the lateral shape of the electron distribution and the height of the 
shower maximum. Due to the statistical nature of shower development, the 
height of the shower maximum is subject to strong fluctuations. Showers, 
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which have started high in the atmosphere show a flat lateral electron distri- 
bution, as electrons in the electromagnetic cascade suffer more from multiple 
scattering processes. Such showers are called old and are characterised by a 
large value of the age parameter. Young showers have started deeper in the at- 
mosphere and had their maximum more close to observation level. This results 
in a steeper lateral electron distribution, which corresponds to a smaller value 
of the age parameter. Apart from fluctuations, the height of the shower maxi- 
mum depends on energy and mass of the shower initiating primary. Therefore, 
the lateral shape parameter is also sensitive to the mass of the primary. 



The mutual interrelation of several independent parameters on which the form 
of the lateral shape depends, makes it a delicate task to draw unique conclu- 
sions from the results of the measurements. Moreover, the interpretation of the 
measured raw data requires a profound understandig on the details of the de- 
tector response functions. Sophisticated simulations of the whole event chain 
that is initiated by the first collision of an ultrahigh energy cosmic ray particle 
with an air nucleus in the upper atmosphere and ends with the registration 
of electronical signals in the various detector components are a prerequisite to 
any reliable analysis of the lateral distributions of all particle components. 



This view encourages to measure the secondary particle components sepa- 
rately, which from an experimental point of view requires several detector 
components to be operated simultaneously. The detector array of the multi 
detector setup KASCADE (Karlsruhe shower core and array detector) [6] is 
designed to disentangle the electromagnetic, the hadronic, and the muon com- 
ponent of the shower disc. Lateral distributions of electrons, hadrons, and 
muons (for different muon threshold energies) in the primary energy range 
5 • lO^^eV < E < lO^'^'eV as measured with KASCADE have already presented 
in a previous paper [7]. In this paper, the lateral distributions of electrons and 
muons in EAS events as measured with the KASCADE array detectors will 
be compared with the predictions of detailed Monte Carlo calculations, which 
comprise the simulation of the full cascading process of EAS, the simulation 
of the array detector response and the final data reconstruction mechanisms. 
Whereas in [7] the parameterisations of the lateral distributions were anal- 
ysed for mean values only, here the reconstruction is also performed on a 
single event basis. Special emphasis is given to investigations of the shape of 
the lateral distributions, the so-called 'lateral age', and its dependencies on 
primary energy and mass of the cosmic rays. Contrary to [7], the hadronic 
component, measured with the KASCADE central hadron calorimeter, will 
not be considered in the present analysis, as well as measurements from the 
additional KASCADE muon devices. 
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The paper is organised as follows: After a brief description of the experimental 
setup an overview on the simulation methods is given. Then we shortly outline 
the data reconstruction scheme. A more detailed explanation is given on the 
method we use for the reconstruction of electron numbers and on function 
used to describe the measured and simulated lateral shapes. This is followed 
by the presentation of mean lateral distributions for muons and electrons, as 
measured with the KASCADE array and a comparison with the results from 
the simulations. Then we show the results for the lateral shape of individual 
showers and its dependence on the shower observables electron and muon 
number. The simulations results, the data will be compared with, are mostly 
based on the hadronic interaction model QGSJet [8], but simulations with 
lower statistics based on the SIBYLL model [9] have also been performed. 
Therefore, results based on the SIBYLL model and the differences in the 
predictions of both models are discussed at the end of the paper. 



2 The KASCADE experiment 

The KASCADE experiment is located at the site of the Forschungszentrum 
Karlsruhe at an altitude of 110 m above sea level. A central hadron calorime- 
ter is surrounded by a rectangular array of 252 scintillation detector stations, 
equally spaced by 13 m and covering an area of 200x200 m^. In addition, 
there is a muon tracking detector with an effective area of 128 m^. The ex- 
periment measures the hadronic, muonic and electromagnetic components of 
extensive air showers in the energy range of 5-10^^ eV up to 10^^ eV of the 
primary particles. A detailed description of the experiment can be found in [6]. 

The 252 detector stations of the KASCADE detector array are organized in 
16 electronically independent clusters. Each cluster consists of 16 stations, ex- 
cept the inner four clusters, where one station per cluster had to give way to 
the central detector. The stations of the inner four clusters contain four liquid 
scintillation detectors, each with an area of 0.8 m^ read out by one photo- 
multipher. The stations of the outer clusters contain two such detectors with 
1.6 m^ total area. All photomultiplier signals of a detector station are added, 
and the integrated charge of the signal is recorded, together with the time of 
the earliest detector hit by a shower particle. These detectors are designed to 
measure arrival time and energy deposits of the electromagnetic component 
of the showers and are therefore referred to as e/7-detectors here. 

Additionally, the stations of the 12 outer clusters house 3.2 m^ plastic scintil- 
lation detectors below a shielding of 10 cm lead and 4 cm of iron, which gives 
0.23 GeV threshold for muons. Each detector is read out by 4 photomultiplier 
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devices and in turn yields time and energy deposit information. Again, the 
sum of the multipUcr signals is recorded together with the hit pattern and the 
time of the earliest detector hit. These detectors measure the muon compo- 
nent and are referred to as muon-detectors here. 

The shower observables which are reconstructed with the KASCADE array 
data are core position, shower direction and the lateral distribution of elec- 
trons and muons. Prom this, the shower size, expressed as total number of 
electrons N^. above 3 MeV and a lateral shape parameter will be derived. For 
the muon component only the total number of muons A^^ above 100 MeV can 
be estimated. Due to the low muon densities, a reliable determination of the 
lateral shape parameter is in general not possible for single event analysis. 



3 Monte Carlo simulation 

A reliable interpretation of the data requires a detailed understanding of the 
physics of shower development, as well as a detailed knowledge of the detector 
response. The whole event chain, starting with the primary interaction in the 
upper atmosphere, followed by the cascading of the shower particles through 
the air up to the response of the detectors at KASCADE ground level has 
been simulated carefully. 

The simulation of extensive air showers is performed with the program COR- 
SIKA (version 6.156 and higher) [10]. For the high energy hadronic interac- 
tions, the models QGSJet (version 01) [8] and SIBYLL (version 2.1) [9] are 
used. Hadronic interactions with energies below 200 GeV are treated with 
the FLUKA code [11], and the electromagnetic component is treated with 
the EGS4 package [12]. The showers were simulated in the energy range from 
1 • 10^^ eV to 1 • 10^'' eV. With respect to computing time the distribution 
in energy was chosen to follow a power law with spectral index ^ — —2. To 
represent different primary masses, the set contains equal numbers of showers 
for five different primary types, namely protons, helium, carbon, silicon and 
iron. The position of the shower cores were distributed randomly all over the 
array. For the shower directions an isotropic distribution was chosen. Output 
of the program is a list of all particles reaching KASCADE ground level to- 
gether with their coordinates, arrival time, 3-momenta and particle type. 

These data are input to the KASCADE simulation program, which is based 
on the GEANT3 package [13]. The simulation covers the whole experiment, 
with all its detectors modelled in great detail. All particles from a CORSIKA 
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simulated shower are tracked through the detectors, the surrounding air and 
the absorber materials. Secondaries, created in interactions with the detector 
materials are likewise followed. In case of the array stations, energy deposits 
and timing information are gathered during the tracking step and converted 
into a photomultiplier signal and a signal time taking into account the light 
coHecting geometry of the detector as well as the specific properties of the 
scintillation material. Output of this program, concerning the array part, are 
arrival times and multiplier signals for e/7- and muon-detectors in exactly the 
same format, that is written by the real experiment after the calibration pro- 
cedure, i.e. the simulated data can be analysed with the standard KASCADE 
reconstruction software. A more detailed description of the KASCADE array 
simulation can be found in [14]. 



4 Reconstruction of particle density distributions and shower pa- 
rameters 

Detector simulations applied to CORSIKA shower events have also been exten- 
sively used for the development and testing of the array data reconstruction 
algorithms and procedures. As the KASCADE array reconstruction scheme 
has already been described in several previous papers [6,7,14], only a brief 
overview will be given in this chapter. Those parts of the analysis chain, which 
are concerned with the reconstruction of the lateral electron distribution and 
its properties will be described in more detail, as they were subject to modi- 
fications applied for the present analysis. 

4.1 Reconstruction of the e/7 component 

Shower direction and shower core position, as well as shower size and the 
lateral form parameter (usually known as age parameter) are reconstructed 
from energy deposits and detector response times of the e/7-detectors using 
an iterative procedure involving three steps. 

In the first step a rough estimate of the shower direction, core position and 
shower size is obtained using fast and robust algorithms, which don't rely on 
any fit procedure. In the second step, the shower direction is determined more 
accurately, by evaluating the arrival times of the first particle in each detec- 
tor. This yields an inclination resolution better than 0.3 degrees for showers 
with IgNg > 4.5 [6]. Then, corrections depending on core position and shower 
inclination are applied to the individual detector energy deposits. From this, 
particle numbers and corresponding particle densities are calculated for each 
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detector. A 4-parameter fit to the spatial distribution of tlie particle densities 
yields core position and shower size, and in addition, a lateral shape parame- 
ter of the charged particle density distribution. The core position resolution at 
this level of reconstruction is better than 0.3 m for showers with Ig A^e > 4.5, 
the shower size resolution at the percentage level [6]. The e/7-detector sig- 
nals contain contributions also from the muon component, for which niTist 
be corrected for, which is performed in the third step of the reconstruction 
scheme: As the analysis of the muon-detector data proceeds in parallel to the 
analysis of the e/7-detector data, the total muon number A^^ is known at the 
time of step three from the step two muon analysis. Therefore, the expected 
muon density can be estimated individually for each c/7-dctcctor. The result- 
ing signal contribution is then accounted for in the detector probabilities for 
the Likelihood minimisation function and the combined lateral density distri- 
butions are fitted. Since the accuracy in core position is in general not further 
improved in this step, only shower size and lateral form parameter are var- 
ied, yielding the final values of the total electron number N^, and the shape 
parameter of the lateral electron density distribution. 



4.. 1.1 Reconstruction of particle densities 

Reconstruction of the lateral particle density distribution requires to interprete 
the measured energy deposits in terms of particle numbers. After correcting 
for different track lengths in the scintillator due to shower inclination, special 
attention must be paid to the 7-component. The electrons are accompanied by 
a multitude of 7-particles, which fake additional electrons, because the photon 
efficiency of the scintillation detectors is roughly 10%. The percentage of fake 
electrons strongly depends on core distance, because the mean 7/e-ratio is 
a function of core distance. In addition, the detector efficiency for electrons 
decreases with increasing core distance, as their kinetic energy distribution 
becomes more and more soft with growing core distance. Moreover, both ef- 
fects depend on shower size and primary particle type. 

A shower size dependent lateral energy correction function (LECF) has been 
derived using the Monte Carlo simulations described in chapter 3. This func- 
tion gives the average expected energy deposit per shower electron as a func- 
tion of core distance and shower size and thereby accounts for the additional 
deposit due to accompanying photons and the dependence of the detector 
efficiency on the kinetic energy of the electrons. Dividing the measured en- 
ergy deposit by the expected deposit per shower electron, an estimate on the 
number of electrons, hitting the detector is given. 

Figure 1 left shows LECFs for proton induced showers for different shower 
sizes. The figure also shows the parametrisation of the mean 7/e-ratio g-y/g 



7 




^ dep 



q t 

50 



JO" 



^ 2 

s 

e 

& 10 



100 150 200 250 300 350 400 
core distance r [m] 



10' 



10" 




- true density 

reconstructed density 



20 40 60 80 100120140160180 200 220 
core distance r [m] 



Fig. 1. Left: Proton LECF's for IgA^g = 5 (dashed lines), IgA^e = 6 (full lines) and 
IgiVg = 7 (dashed dotted lines). Also shown are corresponding 7/e-ratios and mean 
energy deposits for electrons and photons as parametrised and used for the LECF. 
Right: Mean electron densities as reconstructed from simulated detector data for 
proton induced CORSIKA showers. The solid lines give the original CORSIKA elec- 
tron densities. Only electrons and photons have been tracked through the detectors 
in this case. 

and the mean energy deposits E^^^, Ej^^ per electron and photon, respectively. 
The energy E^^ deposited in average by Ug electrons and accompanying 
photons in the detector is then given by 



rptcft 
^dep 



^dep^ 



(1) 



From this the functional form of the LECF ist derived according to 

Jlecf = — ^ = E^^p + • -E'dep ~ -^dep + li/e ■ Ej^p (2) 



where the mean electron and photon energy deposits are assumed to depend 
only on the mean kinetic particle energies. 



The sudden fall of the LECF starting at the shower center reflects the decrease 
of the mean kinetic energy of electrons and 7's with increasing core distance. 
At a core distance of about 30 m, the corresponding loss in detector efficiency 
becomes compensated by the rising 7-electron ratio, which yields an increasing 
fraction of fake deposit due to 7-particles. At large core distances, 7's fake 
nearly half one of the average energy deposit per electron. 

Comparing LECFs calculated for proton and iron primaries, one finds them to 
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differ by less tlian 1% witliin 100 m from the shower core. For larger distances, 
the differences do not exceed 5%. So one is free to use a common LECF for the 
analysis of data, where the primary particle type is unknown. Moreover, the 
variation with shower size is less than 100 m. There is a neghgible dependence 
also on the inchnation of the showers, as geometrical effects are corrected for 
in the reconstruction procedures. 



The right part of Figure 1 shows results of the reconstruction when applied to 
simulated detector deposits, calculated from CORSIKA showers with the de- 
tector Monte Carlo, as described above. To check the reliability of the LECF, 
only electrons and 7-particles have been tracked in this case. The reconstructed 
electron distributions compare well to the distributions of the CORSIKA elec- 
trons. Deviations from the original distribution are found only for large show- 
ers and small core distance. This however is due to saturation effects in the 
e/7-detectors, which are included in the detector simulation. 



4- 1-2 Reconstruction of shower parameters size and lateral shape 

A theoretically motivated function for the description of the lateral electron 
density distribution p(r) is given by the so-called Nishimura-Kamata-Greisen 
function (NKG) [3,4,5] 

— ) (1 + —) . (3) 



with the age parameter s, the Moliere radius tm and the normalising factor 
c{s) 

27rri,r(s)r(4.5-2s)' 



This function has been derived analytically for the case of purely e/7-induced 
air showers but is also used to describe the lateral electron distribution of 
hadron induced showers. It is however known [15,16,17], that the NKG func- 
tions have shortcomings in fitting measured EAS electron distribution, most 
obvious at large core distances. This deficit is usually addressed to the fact, 
that the NKG-function was derived for electromagnetic cascades, whereas 
hadron induced air showers are a superposition of a large number of inde- 
pendent electromagnetic showers. 



In a typical KASCADE event, the detector distances to the core may extend 
up to 200 m. In case of large showers with sizes well above Ig A^e — 6, which 
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Fig. 2. Left: Comparison of fit results for the standard NKG function and the mod- 
ified NKG function. Shown are, as an example, fitted CORSIKA electron distribu- 
tions (i.e. without detector simulation) for proton and iron showers and two different 
shower sizes each. Right: Residuals of the fitted function and the simulated COR- 
SIKA electron distributions for proton and iron induced showers and two different 
shower sizes (dark symbols represent the small, grey symbols the large iVg-bin). The 
iron symbols mostly overlap. 

roughly correspond to a primary energy of 10 PcV, the detectors close to the 
shower core become saturated and must be rejected for the analysis. Thus, the 
lateral fit range differs for small and for large showers significantly in both, 
the upper and lower bound. The deviation of the NKG-function from the true 
shape of the lateral electron distribution therefore gives rise to systematic er- 
rors in the lateral shape parameter (age s in the NKG-formula) , which depend 
on the lateral fit range, and thereby on shower size. 



It has been pointed out [7] , that fixing the age parameter s and instead varying 
the scale parameter Tq (Moliere radius Tm in the NKG-formula) considerably 
improves the fit behaviour of the NKG function. Unfortunately, this does work 
well only for mean lateral distributions constructed from a large number of 
showers. When fitting individual events, which suffer from large statistical and 
physical fluctuations between the detector stations, this method has proven 
to be significantly unstable. Especially for small showers it yields unreliable 
results. 



Other ways to cure the defects of the NKG function for describing electron 
lateral distributions of hadron induced air showers, is to replace it by a different 
function (e.g. as in [18]) or to modify its functional form by changing the values 
of the exponents. Indeed, this gives a better adaption to the shape of the lateral 
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density distribution. For this, we replace equation (2) by 

/ \ s-a / ^ \ s-f3 

P-N,-c(s)-{-) (l + -) , (5) 



with 



5is) = ^-^fl^ (6) 



Testing this function with Monte Carlo data, as optimum values for the ex- 
ponents we have found a = 1.5 and (3 = 3.6 , when Tq = 40 m for the scale 
parameter is used. Optimum in this case means an almost negligible system- 
atic uncertainty in the reconstructed shower size A^e over the full KASCADE 
range as shown in Figure 3. With these values of a,l3 and Tq, equation (5) 
limits the new parameter s to the range —0.5 < s < 1.55. At the same time, 
of course, it looses its numerical relation to the longitudinal development of 
the electromagnetic cascade, as it is often mentioned for its original form. In 
the following, the new parameter s will be called the shape or form parameter 
of the lateral density distribution. 



As an example, the left part of Figure 2 compares both variants of the fit 
function when applied to mean lateral electron distributions derived directly 
(i.e. without detector simulation) from CORSIKA simulated proton and iron 
showers for two different shower sizes. The modified NKG-function adapts 
much better to the lateral shapes, which becomes most obvious at large core 
distances. This can be seen even better at the right part of Figure 2, which 
compares relative deviations of the fit function from the distributions for the 
original NKG-function and the modified one. Apart from the very vicinity 
of the shower core, the modified NKG-function describes the shape of the 
mean lateral distributions with significantly smaller residuals over the whole 
fit range. 



The benefits of the modified NKG-function when applied to individual showers 
are shown in Figure 3. The left part shows results for the estimate of the total 
shower size, derived by both functions. With the modified function, the sys- 
tematic uncertainty in shower size is almost zero over the range 5 < Ig A^g < 7 
and still below 5% between 4.5 < IgA^e < 7.5. Even more convincing, this 
does not depend on the primary particle type, contrary to the results from 
the original NKG function, which fits iron induced shower profiles with a larger 
systematic uncertainty than proton induced ones. 
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Fig. 3. Left: Deviation of reconstructed from true shower size for NKG and modified 
NKG function, when fitting individual showers. Right: Reconstructed shape param- 
eter as a function of shower size using NKG and modified NKG function. The rise 
of the age parameter for showers with Ig iVg > 6.5 is an artefact of the shortcomings 
of the NKG function in describing the lateral distributions at large core distance. 
The vertical dashed lines shows the KASCADE threshold. 

The increasing systematic error towards small shower sizes below IgA'^e = 5 
in the case of iron primaries is related to the strongly growing /i/e-ratio. For 
Ig iVg = 4, the muon density becomes comparable to the electron density even 
at small core distances. This makes it difficult to disentangle both components 
in the e/7-detector analysis and finally leads to an overestimate in shower size. 
For IgA^e > 7 on the other hand, more and more detectors in the vicinity of 
the shower core get saturated. This reduces the available lateral fit range and 
results in a quickly growing underestimate of the shower size. Both effects 
ultimately establish the limits in the primary energy range for this analysis to 
the region between 5 ■ 10^^ eV and 10^^ eV. 



The right part of Figure 3 compares the results of the original and the modi- 
fied NKG function for the lateral shape parameter as a function of the recon- 
structed shower size. As already mentioned, the absolute value of the shape 
parameter is related to the choice of the scale radius rg and the values taken 
for the exponents a and (3 of the modified NKG-function and is therefore 
shifted to smaller values. Apart from that, the results of the conventional fit 
function exhibit a rise of the shape parameter for IgA^e > 6.6 due to the de- 
scribed shortcomings of the original NKG-function and the shrinking of the 
lateral fit range with growing shower size, when detectors near to the shower 
core become saturated. This artefact is absent with the new fit function. For 
shower sizes below the KASCADE threshold the age distributions are effected 
by the selection procedures. 
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Fig. 4. Muon lateral densities as reconstructed from simulated muon-detector data 
compared to the true CORSIKA densities. 

4.-2 Reconstruction of the ii-component 



The analysis of the muon-detector signals follows the same line as in the case 
of the e/7-detectors. Again, the raw signals are subject to corrections due to 
shower inclination and core position using a corresponding muon LECF, which 
has a much simpler form than in the e/7-case. Additionally, corrections for 
punch through from the e/7- and hadron component are applied before the 
muon densities for every detector are calculated. Detectors closer than 40 m 
to the shower core must be excluded from the analysis, because in this region 
punch through dominates the signal. 



The total muon number is estimated by fitting a modified NKG function with 
exponents a = 1.5 and f3 = 3.7. For the muon component this gives only a 
moderate improve over the original NKG-function, which is known to fit muon 
lateral distributions already quite well, provided the scale parameter is chosen 
appropriately. Here we take Tq^ = 420 m. Due to the low muon densities, a 
2-parameter fit on a single air shower basis proves unreliable. Therefore the 
total muon number Nfj_ of the shower is estimated with a fixed lateral form pa- 
rameter Sn and a 1-parameter fit. The muon lateral shape is parametrized 
as a function of shower size N^. from CORSIKA simulations and the actual 
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value is chosen event by event during the iterative reconstruction. For the con- 
sidered data sample varies between 0.81 and 0.75, slowly decreasing with 
increasing shower size. 

Reconstruction results for muon densities are illustrated in Figure 4, where 
muon distributions are compared for simulated showers with the corresponding 
true CORSIKA distributions for several ranges of the reconstructed total muon 
number N^. For individual showers, the accuracy of this observable is moderate 
compared with shower size A^e and is not better than 10 to 20%. This results 
from the poor muon statistics in case of small showers and punch through at 
large shower sizes. The muon size iV^ is input to the correction of shower size 
and shape parameter in the third step of the e/7-detector analysis. 



5 Compcirisons of KASCADE data with Monte Carlo simulations 

For the results presented here all measured showers with log A'e > 6 have been 
taken into account. This sample comprises about 170 000 events measured in 
a period of nearly 8 years. Additionally, these data have been enriched with a 
sample of the many small showers, which hit the array more frequently due to 
the steep energy spectrum of cosmic rays. For this, two KASCADE runs with 
in total about 2.5 million recorded shower events were added to the data set. 

All showers, real or simulated, included in this analysis were subject to the 
same reconstruction procedure and to the same cuts concerning trigger condi- 
tion, core position, inclination angle and shower shape parameter. Showers are 
restricted to core distances less than 90 m from the array center and zenith 
angles less than 30 degrees. An additional cut for showers with shape parame- 
ter value s > 1.4, which is close to the upper boundary of the mathematically 
possible range, excludes showers which are frequently misreconstructed inside 
the array but actually had its core outside, or are just very small showers 
which fluctuated in such a way that the reconstruction overestimated their 
size by a large amount. Indeed, showers of the second kind are already very 
efficiently excluded, by comparing shower sizes as reconstructed during step 
one and step two of the analysis, and cutting on those events, where the differ- 
ence in the estimated sizes exceeds the expectations due to the uncertainties 
in both methods considerably. Since the energy distribution of the simulated 
showers represents a spectral index of 7 = —2 while the real data follow an 
index with —2.6 < 7 < —3.1 fluctuations to larger values in shower size would 
be less pronounced in simulations. Therefore appropriate statistical weights 
have been given to the simulated events, as will be explained below. 
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The analysed KASCADE data set is first compared to the predictions of tlie 
QGSJet model, based on a sample of 1.7 Million simulated events. In addition 
a set of showers with half that statistics but based on the SIBYLL model was 
analysed (section 5.4). 



5. 1 Lateral distributions of muons 

Figure 5 shows mean lateral distributions from the simulations based on the 
QGSJet model and compares them with the data. Showers have been sampled 
in ranges of the reconstructed total muon number per single event. Because of 
the steeper energy spectrum of the data, fluctuations in muon number would 
give significantly larger contributions to higher iV^-bins than they would for 
simulated events. To account for this effect and also to show its result on the 
form of the lateral distributions, the simulated events have been analysed for a 
statistical weight representing a spectral index 7 = —2.6 as well as for weights 
giving an index 7 = —3.2. The spectral index of the data varies with energy, 
but will lie somewhere between these values. The small shaded bands in Figure 
5 give the simulation results within this bounds. The lower bound of the band 
always corresponds to the larger absolute index value, i.e. to 7 = —3.2. The 
width of the band for an individual primary mass is in the order of the symbol 
size, only. It is obvious, that the form of the lateral distribution function 
does only weakly depend on primary energy. Moreover, it can be seen, that 
the bands of proton and iron nuclei at least partially overlap. The simulations 
therefore predict, that the form of the muon distribution is not sensitive to the 
nature of the primary particle. At low energies, proton induced showers show 
slightly steeper lateral shapes compared to iron primaries, but these differences 
vanish for higher energies. Even for small showers the differences in the density 
distributions for showers of either type do not exceed ten percent. 

The lateral distributions derived from the data agree quite well with the Monte 
Carlo predictions. The simulated distributions describe the measurements over 
the full KASCADE range of core distances and primary energies. The figure 
thereby shows clearly, why the muon measurements at KASCADE are sensi- 
tive to the primary energy, but give no valuable information on the elemental 
composition of cosmic rays. This is also found by using the SIBYLL model. 



5.2 Lateral distributions of charged particles 

Figure 6 presents lateral distributions of charged particles as measured by the 
e/7-detectors together with distributions derived from simulations. Though 
the bulk of particles are electrons, also muons contribute to the e/7-detector 
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Fig. 5. Lateral density distributions of muons as measured with the KASCADE 
array muon-detectors. The shaded bands cover the range of the Monte Carlo simula- 
tion results with respect to five different elemental masses including an uncertainty 
in the spectral index within the range —3.2 < 7 < —2.6. For reasons of clarity, only 
the results for two elemental masses are shown additionally. For these, a spectral 
index 7 = — 3 is assumed. 

signals, while contributions of photons are corrected by the LECF, and con- 
taminations by hadrons are negligible. Therefore, the distributions presented 
here contain besides electrons also muons. All showers have been grouped 
according to their reconstructed total muon number. 

The shaded bands of Figure 6 again indicate the range of uncertainty, which 
results, when the spectral index is varied between —3.2 < 7 < —2.6 and the 
primary mass from proton to iron. Again, the lower bounds corresponds to the 
larger absolute value 7 = —3.2, the upper bounds to 7 = —2.6. In addition, 
within each band the lateral distributions for five different primary masses 
are shown, assuming a spectral index of 7 = — 3 for each. The curve close to 
the upper bound of each band now results from the lightest primary, the one 
close to the lower bound belongs to iron. This expresses the well known fact 
[19], that showers originating from light primaries are more electron rich at 
sea level than showers induced by heavy nuclei. The e/yU-ratio is therefore a 
common starting point for the analysis of the chemical composition of cosmic 
rays. Furthermore, the figure illustrates, that showers induced by light pri- 
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Fig. 6. Lateral density distributions of electrons including muons as measured with 
the KASCADE e/7-detectors and by simulations. The shaded bands cover the range 
of the simulation results with respect to five different elemental masses including 
an uncertainty in the spectral index within the range —3.2 < 7 < —2.6. The results 
for each elemental mass, which are also shown, assume a spectral index 7 = —3. 

maries are predicted to exhibit a slightly steeper electron distribution than 
showers stemming from heavy nuclei. 



Comparing real data with simulations, one sees that small showers with Ig A^^ 
4 fit simulated proton and helium distributions quite well, while large showers 
with IgiV^ pa 7 arc best described by the silicon and iron distributions, i.e. 
by heavier primary particles. The figure visualizes the known [6,20] variation 
in the average e///-ratio with shower size, which indicates, within the scope 
of the QGSJct model, a transition of the primary particle composition from 
light elements at energies below to heavy elements at energies above the knee. 
A detailed analysis of the electron-muon-number frequency spectrum as mea- 
sured by KASCADE is subject of a separate paper, focused on the chemical 
composition of cosmic rays [21]. 



A closer look to Figure 6 also reveals discrepancies between the data and the 
QGSJet Monte Carlo predictions. At low energies, i.e. small muon numbers. 
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the shapes of the measured distributions are shghtly flatter than those of the 
simulated proton and hehum distributions. For the smallest showers consid- 
ered here, the measured particle densities even exceed the expected range of 
densities at core distances beyond 120 m. This kind of shortcoming can not be 
cured by any assumption on the elemental composition and will be discussed 
in more detail in the next chapter. 

Looking at the highest energies, where data are best described by heavy pri- 
mary showers, simulations seem to show slightly lower densities and a flatter 
lateral behaviour at small core distances. In this region however, data suffer 
severely from overflows, from which additional uncertainties result. 

It might be worth mentioning, that the shortcomings described here must orig- 
inate completely from processes involved in the generation and development of 
the e/7-component. The muon component contributes at low shower energies 
considerably to the lateral distributions measured with the e/7-detectors. But 
the form of the muon lateral distribution is well described by the Monte Carlo 
simulations and therefore can not be drawn to explain the observed deviations. 



5.3 The lateral shape parameter 

The most obvious differences between the lateral shapes of the individual el- 
ements shown in Figure 6 are simply the amplitudes of the density functions, 
and are related to different e//x-ratios. A more subtle quantity, which in this 
kind of representation is difficult to compare, concerns the functional form or 
the slope of the lateral distribution. 

The relations of shower sizes and with the shape parameter s and with 
the primary mass are illustrated by Figure 7. It shows the mean reconstructed 
shape parameter value as a function of the observables IgA^e and IgA*"^ for 
QGSJet based simulated showers. The individual events are weighted in en- 
ergy to represent an elemental composition according to the results for the 
QGS modell as described in [21]. The lines overlaid to the distribution repre- 
sent linear approximations to lines of maximum probability for showers of a 
single element (here shown only for proton and iron) but variable energy. It 
is obvious, that showers from light primaries are younger in average, i.e. have 
smaller shape values compared to heavy primaries, and showers of high energy 
are younger than low energy ones. 
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Fig. 7. Mean shower shape parameter s as a function of Ig A'^ and Ig for simulated 
showers based on the QGSJet model with a composition adopted from [21]. Also 
shown are the lines of maximum probability for proton and iron induced showers 
and the coordinate system, used for the comparison of simulation results and data. 

The lines of maximum probability in Figure 7 offer one axis of a natural, 
rectangular chosen coordinate system to compare the data with the simula- 
tion results. The new coordinates are related to the Ig A'^^ — IgiVg— system by 
a simple rotation around the origin with an angle of 51,6 degrees, obtained 
from the simulations. This coordinate system simply adapts the form of the 
event distribution in the IgA^^ — \gNe— plane and will be used in the follow- 
ing to compare data and Monte Carlo results only on the basis of measured 
(or simulated) observables. While the coordinates along the lines of maximum 
probability may be associated with an energy estimator £igE, the coordinates 
perpendicular to this direction measure a mass estimator These are surely 
not the best possible estimators for energy and mass but we do not want to 
draw quantitative conclusions from them. Therefore the numerical values of 
the new coordinates are simply the values obtained by the rotation from the 
old IgA^^ — IgiVe— values. 

The shape parameter as a function of the energy estimator e\^E is shown in 
Figure 8 for both, data and Monte Carlo simulations of all five elements. The 
Monte Carlo results confirm that showers induced by heavy primaries are older 
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Fig. 8. Reconstructed shape parameter (of the modified NKG function) as a function 
of the energy estimator eigE for KASCADE data, five primary masses and a model 
composition as simulated using QGSJet. The scale on top gives a rough estimate 
for Ig E in GeV. 

compared to showers of light primaries. With increasing energy the shape pa- 
rameter value decreases for all simulated elements and reflects the fact, that 
the height of the shower maximum decreases with increasing energy. 



The data fit into this picture only qualitatively. Up to an energy of about 
10 PeV, they follow the line of carbon. For low energies, this suggests a rela- 
tively heavy composition, which clearly disagrees with the predictions of figure 
6. For higher energies, the lateral shape parameter stays almost constant and 
crosses the hne of iron at an energy of about 30 PeV. Beyond this crossing 
point, the absolute values of the measured shape parameter cannot be ex- 
plained by any elemental composition within this Monte Carlo model. 



For a more detailed investigation the data distributions are compared with 
what would be expected from the simulations, once a reasonable elemental 
composition is given. For this, the simulated shower events of the five ele- 
mental masses have been weighted with individual energy spectra, which have 
been reconstructed from an analysis of the measured A^eZ-^^i-spectrum using 
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Fig. 9. Reconstructed shape parameter (of the modified NKG function) as a func- 
tion of the mass estimator e„j for K AS CADE data and five primary masses as 
simulated using QGSJet for different ranges of the energy estimator eig^;. The fine 
of maximum probability for protons corresponds to a value Em = —0.22 and for iron 
it is Em = —0.45. 

a sophisticated unfolding algorithm based on the same model QGSJet. The 
resulting model composition favours light elements before the knee and a sig- 
nificant contribution from heavy elements at energies above the knee [21]. 
The effect on the shape parameter as a function of the energy estimator is 
also shown in Figure 8. It is remarkable, that the hne of the measured shape 
parameter values runs almost parallel to the line representing these adapted 
Monte Carlo predictions, but is displaced by a nearly constant amount of 
As ~ 0.05 over the whole energy range. 



The behaviour shown in Figure 8 may therefore be interpreted in the same 
way as the form of the lateral shapes in Figure 6. The almost constant value 
of the shape parameter for energies beyond 10 PeV can be understood as the 
result of a transition from light to heavy nuclei in the elemental composition 
of cosmic rays. The offset between the lines of measured and simulated shape 
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simply states, that the simulations in general yield slightly steeper shapes than 
observed in real showers. 

This can be seen even more clearly when looking along the lines of constant 
values of the mass estimator This view is given in Figure 9 for several 
ranges of the energy estimator Ei^e, i-e. slices of Figure 7 perpendicular to the 
lines of maximum probability. Here higher mass values correspond to smaller 
shower size and larger muon content, i.e. to showers which are electron-poorer. 
Remarkably, all elemental masses can be seen to follow the same (energy de- 
pendent) functional dependence between the shape parameter and the mass 
estimator. On the other hand this may be expected, because showers of heavy 
primaries develop higher in the atmosphere compared to showers of light nuclei 
but same energy. However, the total number of electrons present in the shower 
maximum is, for showers of the same energy, roughly independent of the kind 
of the primary nucleus. Therefore, a deeply penetrating iron shower may not 
be distinguishable in shape from a proton shower, which developed very early 
in the atmosphere. The shape of real showers however does not follow this 
functional dependence. Measured showers are older in average, with values 
deviating by an amount As ~ 0.05 from simulated showers, with increas- 
ing tendency at the highest energies. This might indicate, that real showers 
develop at higher altitudes than predicted by the simulations or/and that 
multiple coulomb scattering plays a distinct and more pronounced role in the 
development of the real electromagnetic cascade expected from simulations. 



5.4 Comparison with the SIBYLL model 

The data have also been compared with simulations based on the SIBYLL 
model. While the e//x-ratio for SIBYLL showers is larger in general, no no- 
table differences to the QGSJet model were found when comparing the shapes 
of the lateral distributions of muons for equal total muon numbers. The lateral 
shapes of the electron component show only small differences compared to the 
QGSJct model, as can be seen from Figures 10 and 6. The SIBYLL calculated 
shapes predict a more heavy composition, as a result of small differences in 
the e///-ratios. In addition, the mean lateral electron distributions appear a 
bit younger. 

Investigating the dependence of the lateral form parameter on the energy esti- 
mator as done in section 5.3, one can see in Figure 11 that SIBYLL describes 

the data worse compared to QGSJet. The SIBYLL iron curve crosses the 
data already at an energy of about 10 PeV, so there is no explanation for 
the measured shape values within this model for larger energies. Comparing 
with Figure 8 one finds that the mean shape of SIBYLL showers in general is 
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Fig. 10. Same as Figure 6, but with SIBYLL generated showers, 
smaller by As ~ 0.05 compared to QGSJet. 

It may appear surprising then, that individual SIBYLL showers follow the 
same functional dependence on the mass and energy estimators (and there- 
fore also on IgA'^e and IgA*"^) than QGSJet showers do, as can be seen from 
Figure 12. This shows that the longitudinal development of the electromag- 
netic component must be very similar in both models. The difference in the 
mean lateral shape parameter is therefore simply due to a different distribu- 
tion of SIBYLL events in the — A^^u— plane (see also [21]). SIBYLL showers 
are more electron rich and produce less muons. The abundancy maximum for 
a given primary energy is therefore shifted to lighter mass values, enhancing 
the weight of younger showers when averaging on the shape parameter. 



6 Summciry and conclusions 



Lateral electron and muon density distributions of air showers as measured 
with the KASCADE array have been compared to the results of Monte Carlo 
simulations based on the CORSIKA program using EGS4 and the two hadronic 
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Fig. 11. Same as Figure 8, but with SIBYLL generated showers, 
interaction models QGSJet and SIBYLL. 



Muon lateral distributions measured with the KASCADE array muon detec- 
tors were found to be well described by the Monte Carlo simulations and 
no significant differences were observed between the two hadronic interaction 
models QGSJet and SIBYLL. Moreover, muon lateral distributions appear 
very similar in shape, independent of the nature of the primary particle, so 
that details of the chemical composition can not show up in the comparison 
of data and simulations. 



Deviations from the Monte Carlo predictions are found for the lateral distri- 
butions of charged particles, which were reconstructed from the measurements 
of the KASCADE e/7-detector array. Common to both models is that they 
suggest a transition from light to heavy nuclei in the chemical composition of 
cosmic rays in the energy range of 1 PeV to 100 PeV. This is consistent with 
the results of an independent study [21] based on a detailed analysis of the 
Ig A^e — lg-^M~ frequency spectrum of KASCADE events. 



Investigating in detail the shape of the lateral distributions, the absolute val- 
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Fig. 12. Reconstructed shape parameter (of the modified NKG function) as a func- 
tion of the mass estimator Em for QGSJet and SIBYLL based simulations and dif- 
ferent ranges of the energy estimator £\gE- For reasons of clarity only two elemental 
masses are shown. 

ues of the measured shape parameter however were found to disagree with 
the predictions of either of the two hadronic interaction models. While both 
models yield the same dependence of the average shower shape on Ig A^'e and 
IgA^"^, the absolute values appear smaller than the measured shape values for 
the whole considered range in IgA^e and IgA*"^. 



Looking at the mean shape parameter as a function of primary energy, SIBYLL 
yields smaller values compared to QGSJet. The reason for this difference 
between the two models is a kind of "lighter" distribution of events in the 
Ig A^e — Ig-^/i— plane in case of SIBYLL: showers of same primary type and 
energy contain on average less muons but more electrons and therefore make 
up a smaller average value for the shape parameter compared to QGSJet. 
However, also the QGSJet predictions still underestimate the results from 
measurements by an almost constant amount of As ~ 0.05. 
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Summarizing, both models are not able to describe the measured lateral dis- 
tribution of the e/7-component correctly. The details of the form of the lateral 
distribution depend on the hadronic interaction mechanism as well as on elec- 
tromagnetic cascading processes. Thus, a variant of the QGSJet model that 
predicts a larger e///-ratio, would give better consistency with data. However, 
the discrepancies might also be burried in the electromagnetic cascading al- 
gorithm EGS4 and its treatment of the multiple coulomb scattering process. 
Further improvements in the simulation models may be necessary, to under- 
stand and remove the discrepancies between data and simulations. Meanwhile, 
these results hopefully may help to stimulate this process and provide some 
additional hints. 
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